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Abstract 


A higher-order bending theory is derived for laminated composite and sandwich beams. 
The recent {1,2} -order theory is extended to include higher-order axial effects without 
introducing additional kinematic variables. This is accomplished by assuming a special 
form for the axial and transverse displacement expansions. An independent expansion is 
also assumed for the transverse normal stress, improving the transverse strain and stress 
predictive capability of the theory. Appropriate shear correction factors based on energy 
considerations are used to adjust the shear stiffness, thereby improving the transverse 
displacement response. A set of transverse normal correction factors is introduced, 
leading to significant improvements in the transverse normal strain and stress for 
laminated composite and sandwich beams. A closed-form solution to the cylindrical 
bending problem is derived, demonstrating excellent correlation to the corresponding exact 
elasticity solutions for a wide range of beam aspect ratios and commonly used material 
systems. Accurate shear stresses for a wide range of laminates, including the challenging 
unsymmetric composite and sandwich laminates, are obtained using an original corrected 
integration scheme. For application of the theory to a wider range of problems, guidelines 
for finite element approximations are presented. 
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1. Introduction 


The high performance aircraft and spacecraft currently being considered for the next 
generation of aerospace vehicles require the use of strong yet lightweight materials. 
Laminated composites, such as graphite fiber reinforced epoxy, are commonly used for 
spacecraft and aircraft where weight is critical. While composite laminates provide higher 
stiffness, strength, and reduce weight over conventional metallic structures, sandwich 
structures can further reduce the structural weight and improve thermal performance 
without sacrificing stiffness and strength. A sandwich structure typically consists of two 
laminated composite or metallic face sheets with a core of foam, metallic honeycomb, or 
other lightweight material. The face sheets provide axial stiffness while the core acts as a 
shear web to carry the transverse shear load. Core materials are typically lighter and 
several orders of magnitude more compliant than the face sheet materials. This presents 
an analytical challenge with respect to the prediction of strain and stress distributions 
through the thickness. The knowledge of accurate, detailed strains and stresses is critical 
to the design of lightweight aerospace structures. Without accurate strain and stress 
predictions, higher factors of safety are often used, thereby increasing the weight and cost 
while reducing performance. 

Numerous bending theories to analyze the response of beam, plate, and shell structures 
have been proposed. Many of the significant developments in elastic theory are 
summarized in review papers by Reissner (1985), Reddy (1989), and Noor and Burton 
(1989). The following discussion reviews the most pertinent developments in beam 
theory, and also makes comparisons to similar plate and shell theories. 

Classical Bemoulli-Euler beam theory is the earliest and simplest approximation used 
for analysis of homogeneous beams. It dates back to 1705 and precedes the theory of 
elasticity by over 100 years (Love, 1952). In the classical bending theory of beams, the 
beam cross section is assumed to be much smaller than the length of the beam and the 
displacements are small compared to the thickness of the beam. The transverse normal 
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and shear strains are assumed to be negligible. As a result, the cross section to the 
midplane is assumed to maintain a constant thickness and remain straight and normal to 
the midplane after deformation. While this approximation is appropriate for thin, 
homogenous beams, neglecting transverse shear deformation causes the deflection of the 
beam to be underpredicted for thick beams. To remedy this deficiency, several theories 
have been developed which account for transverse shear deformation. 

The first beam theory to include transverse shear deformation was proposed by 
Timoshenko (1921). The cross section is no longer restricted to remain normal to the 
midplane. The cross section is assumed to remain plane as in classical theory, but 
transverse shearing is permitted to occur. By accounting for shear deformation, the 
theory provides more accurate response predictions for thin and moderately thick 
homogeneous beams. Reissner (1945) and Mindlin (1951) later extended Timoshenko 
theory, also known as first-order shear deformation theory, to the analysis of 
homogeneous elastic plates. Since then, numerous displacement-based, stress-based, and 
mixed formulation bending theories have been developed, with more recent emphasis on 
the analysis of laminated composite and sandwich structures. 

Although classical and Timoshenko beam theories were originally developed for 
homogenous beams, they can also be extended to the analysis of heterogeneous beams. 
This is accomplished by modeling the laminate using either (1) a layer- wise theory or (2) 
a single layer theory. Both theories assume the laminate to be an assembly of 
homogeneous, anisotropic plies perfectly bonded at the ply interfaces. For layer-wise 
theories, the displacement components are assumed to be piece-wise smooth through the 
thickness. This implies that the function is continuous through the thickness of the beam, 
but at the ply interfaces, the slope of the function may not be continuous. These theories 
are especially useful for the analysis of moderately thick to thick laminates with span to 
thickness ratios ranging from 10 to 4. In single layer theories, the laminate is treated as an 
equivalent single layer, where the displacement assumptions represent some weighted 
average distribution through the thickness. Although layer-wise theories are capable of 
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accurate modeling of laminated composites and sandwiches, the number of independent 
variables is directly proportional to the number of plies in the laminate. Solutions to 
these equations can become computationally expensive and intractable as the number of 
layers in the laminate increases (Reddy and Liu, 1987). The single layer theories are found 
to be the most computationally efficient and will provide the background for the present 
higher-order theory. 

For single layer, displacement-based theories, the displacement approximations are 
assumed to take on a particular polynomial form through the laminate thickness. Each 
term in the expansion adds an additional power of the thickness coordinate, while the 
expansion coefficients are functions of the x coordinate only. These approximations are 
the basis for the development of the strain and stress quantities, and hence will govern the 
complexity and accuracy of the theory. The axial displacements are expanded with a 
polynomial of degree, m, while the transverse displacement may be of a different degree, 
n. Thus, the notation {m,n} may be used to distinguish between the different order single 
layer theories (Tessler and Saether,1991). 

One of the most commonly used approximations for single layer displacement-based 
theories is {1,0} (or first-order) shear deformation theory, i.e., having a linear axial 
displacement and a constant transverse displacement. With the advent of advanced 
composite materials, there has been a great deal of work in the development of plate 
theories which extend the Timoshenko beam theory and Reissner-Mindlin plate theory to 
the analysis of laminated, heterogeneous structures. Stavsky (1959) is apparently the 
first to develop such a theory, and improvements upon this theory have been made by 
Yang, et al. (1966) and Whitney and Pagano (1970). Even with modifications to the 
theory, the detailed stress distributions show little improvement over the classical theory 
for moderately thick heterogeneous beams. Although the transverse shear stresses are 
included in these theories, the transverse normal stresses are neglected. Both transverse 
shear and transverse normal stresses become significant for moderately thick to thick 
laminated composite and sandwich beams. They often contribute to delaminations of the 
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plies and subsequent structural failure. Thick composites and sandwiches generally 
exhibit nonlinear displacement and strain distributions through the thickness which cannot 
be predicted using linear displacement assumptions. Because of the inadequacies of the 
first-order shear deformation theories, higher-order theories have been developed which 
attempt to predict more accurately the response of heterogeneous structures. 

Higher-order theories in this context refer to the class of theories, for either beams, 
plates, or shells, in which the polynomial expansions of the displacement field are of 
higher order than {1,0} approximations. While the {1,0} Timoshenko beam theory 
includes transverse shear deformation, the transverse normal strains are still assumed 
negligible. Essenburg (1975) proposed a { l,2}-order beam theory which includes 
transverse normal effects, and a similar plate theory has been proposed by Whitney and 
Sun (1974). Lo, et al. (1977) developed a {3,2} laminated plate theoiy which produces 
adequate predictions for inplane displacement and stress distributions. The theory has a 
higher degree of complexity having eleven kinematic variables, and its predictive 
capability for the transverse normal and shear stresses is rather poor. 

As a compromise between accuracy and computational efficiency, Reddy (1984) 
developed a {3,0} -order theory. By taking a special form for the inplane displacement 
components, the theory contains the same number of kinematic variables as the first- 
order shear deformation theory. This special displacement form also results in the 
parabolic distribution of the shear strain through the thickness. The shear stresses satisfy 
the traction-free boundary conditions on the surfaces of the plate, an important physical 
condition which previous theories did not enforce. Phan and Reddy (1985) investigated 
solutions to the aforementioned theory, showing improvements over classical theory. 

Tessler and coworkers (1991-1995) developed a { 1,2} -order theory for beam, plate, 
and shell analyses. The theory assumes a special form for the transverse displacement to 
obtain the ‘average’, parabolic shear strain distribution and shear traction free surface 
conditions. The theory is novel in that the transverse strains are not directly derived from 
the strain-displacement relations. The transverse normal and shear strains are assumed as 
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independent polynomial expansions. They are required to be least squares compatible, 
through the laminate thickness, to the strains obtained directly from the strain- 
displacement relations. The resulting thickness distributions for the transverse stresses 
and strains show adequate correlation with results given by elasticity theory, an 
improvement over previous higher-order theories. Tessler (1993a) improved the theory 
further for application to composites by introducing an independent polynomial 
assumption for the transverse normal stress in place of the transverse normal strain. The 
improved theory results in a more accurate representation of transverse normal stresses 
and strains, and is further substantiated by solutions given by Schleicher (1994). The 
{1,2} theory retains the simplicity of the first-order shear deformation theory as far as 
the engineering boundary conditions are concerned. Furthermore, the theory gives rise to 
finite element formulations that are compatible with the first-order shear deformation 
elements. 

Application of the {1,2} theory generally results in excellent predictions for thin and 
moderately thick homogeneous and laminated beams, but the theory has some deficiencies 
in modeling the response of sandwich beams. The linear axial displacement assumption 
cannot model the nonlinear thickness distributions of the axial displacement and strain. In 
thick laminates, this generally results in underestimation of the axial stress, typically the 
largest stress which governs the design of the structure. The other deficiency is the 
transverse normal stress violation of traction conditions on the top and bottom surfaces 
of sandwich beams. 

The main objective of this report is to expand upon Tessler’ s {1,2} theory by 
proposing a higher-order theory of order {3,2} which includes nonlinear axial effects. The 
theory is expected to model accurately the nonlinear axial displacement and strain through 
the thickness of the beam, resulting in a more accurate axial stress prediction. A special 
form for the cubic axial displacement field is used such that no additional kinematic 
variables are introduced into the theory. This enables the present higher-order theory to 
retain the simplicity of the {1,2} theory while improving predictions of the axial response 


5 



of the beam. The transverse normal stress and transverse shear stress predictive 
capability will also be improved for sandwich beams. 

In Section 2, the mathematical formulation of a higher-order theory is presented. The 
{3, 2} -order displacement field is assumed. These displacement approximations account 
for the nonlinear variations of stress and strain quantities typically present in thick 
composite and sandwich structures. As in the {1,2} -order theory, a special form of the 
quadratic transverse displacement is assumed. An independent expansion is also assumed 
for the transverse normal stress, as in Tessler (1993a). However, the shear strain is 
computed directly from the strain-displacement relations. Because of the special form for 
the axial displacement assumption, the traction-free shear stress boundary conditions are 
satisfied exactly. The principle of virtual work is used to derive the equilibrium equations 
and boundary equations for the theory. The hierarchical form of the {3,2} -order 
displacement approximations permits a straightforward reduction to lower-order beam 
theories of Tessler and Timoshenko. 

In Section 3, transverse shear and transverse normal correction factors are determined 
from the energy and traction equilibrium considerations. Furthermore, accurate piece- 
wise smooth shear stresses are determined by integrating the two-dimensional equilibrium 
equation of elasticity theory. A correction procedure is also developed to improve the 
accuracy of this approach for unsymmetric and sandwich laminates. 

In Section 4, an analytic solution to the cylindrical bending problem is presented. 
Numerical results are presented for commonly used aerospace material systems. 
Comparisons are made to the {1,2} theory and three-dimensional elasticity solutions. In 
Section 5, guidelines for finite element formulation and implementation based on the {3,2} 
theory are presented. Finally, Section 6 presents conclusions and recommendations for 
future work. 


6 



2. Mathematical Formulation 


This section discusses the derivation of a higher-order bending theory for 
heterogeneous beams. A displacement field of order {3,2} is assumed and the strains, 
stresses, equilibrium equations, and boundary conditions are derived. 


2.1 Higher-Order Theory for Laminated Beams 

Consider a straight, linearly elastic beam laminated with N orthotropic plies subject to 
the loading shown in Figure 2.1. The beam has a span L and a rectangular cross-section of 
thickness 2h and width b. The orthotropic plies are stacked from the bottom face (z= -h) 
such that the material properties, in general, are functions of the z coordinate. The 
loadings q + and q' are stresses applied normal to the top and bottom faces of the beam 
and may vary in the x coordinate. T i0 and T iL ( i = x,z ) are arbitrary stress components 
prescribed at the ends of the beam. 


z,u 2 






T x0 1 

TT 


■f 

!h 

!|h 



h b H 



Figure 2. 1 Sign Convention for Beam 
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2.1.1 Kinematic Displacement Assumptions 


The first step in developing a displacement-based theory is to introduce kinematic 
assumptions for the displacements. These assumptions are made by first choosing a 
polynomial form for the displacement components. For thick laminated composite and 
sandwich beams, the displacement components are piece-wise smooth and nonlinear 
through the thickness. This contrasts with the predominantly linear displacement 
distributions for thin beams. To capture these higher-order deformation effects, a cubic 
polynomial is assumed for the axial displacement u x and a quadratic distribution is 
assumed for the transverse displacement u z , i.e. 


u x (x, z) = U 0 (x)+ u,(x)? + u 2 (x) C 2 + u 3 (x) ? 3 
u z (x,z) = w(x)+w,(x)^ + w 2 (x) (? 2 +c) 


where ^ — — e [—1,1] is the dimensionless thickness coordinate such that £ = 0 defines 

the midplane of the beam. The four U; coefficients in the axial displacement expression 
are independent unknowns to be determined. The w ; coefficients in the transverse 
displacement are kinematic variables identical to those defined by Tessler (1991). The 
constant C is included in the transverse displacement equation to allow w(x) to represent 
a weighted average transverse displacement. 

Three conventional kinematic variables are introduced and defined, as in Reissner 
(1945), as weighted average quantities through the thickness, such that 
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u ( x ) = ^j_ h u >x x >z)dz 

0 ( X ) = 2M-h U *( X ’ Z > zdz 

w(x) = ^rj. h h u z (x, z)(l-C 2 )dz 


( 2 . 2 ) 


where u(x) is the midplane displacement along the x axis, 0(x) is the rotation of the 
normal about the y axis, and w(x) is the weighted average of the transverse displacement. 
To determine the constant C which satisfies the weighted average displacement 
distribution through the thickness, the transverse displacement expansion from (2.1) is 
substituted into the expression for w(x) in equation (2.2) to yield 

w(x) = irf h [ w ( x ) +w,(x) c+ w 2 (x) (C 2 +c)](l -C 2 ) dz (2.3) 

from which C is found to be equal to ~/ 5 . Similarly, the expressions for u(x) and 0(x) 
from (2.2) are used to determine two of the four coefficients for the axial displacement 
distribution: 


U( x ) = rr J> o( x )+ u i(x)C + u 2 (x)C 2 +u 3 (x)p]dz 
0(x) = fr C l u o( x ) + u.( x ) C + Uj ( x ) C 2 + u 3 ( x )C 3 ] z dz 


(2.4) 


From the above equations, the coefficients u 2 and u 3 are calculated to be 


u 2 (x) = 3 (u(x) - u 0 (x) ) , u 3 (x) = |-(h 6(x) - u, (x) ) 


(2.5) 



The remaining coefficients are determined by using the condition that the shear stress at 
the top and bottom faces of the beam must vanish, i.e. T xz (x,z=±h)= 0 . Because of the 

direct relationship between stress and strain in beam theory, T xz = C 55 y xz , the shear 
strains at the top and bottom faces must also vanish: 

Tx Z ( x > z =± h )=° (2.6) 


Using the strain-displacement relation for the transverse shear strain, y xz = u x z + u z x , 
and equation (2.5) yields 


Y„(x,z) = w(x)+ u - (x) + W ' (x V + 6 ( U(X) ~ 2 U ° (X) h 
> h h 

5(h0(x)-u,(x))z 


+ 


h 3 


•+w 2 (x) J 


U 1 


\h 


2 5 


(2.7) 


where a comma),) denotes partial differentiation. After enforcing the two boundary 
conditions from (2.6), the coefficients u 0 and u, are found to be 


u 0 (x) =u(x) + 


h w ,(x) 


u,(x) 


h (25 0(x) + 5 w(x) x + 4 w 2 (x) x ) 


20 


( 2 . 8 ) 


Substitution of the coefficients u 0 , , u 2 , u 3 , and the constant C into equation (2.1) 

yields the final form of the displacement components: 


u x (x,z) = u+ hCO -£(3 C 2 - 1 ) h w, x -h C (K 2 — 0y 
u z (x,z) = w + t>, +| 2 — 0w 2 
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with the shear angle y defined as 


Y (x) = f [B(x)+ w x (x)]+ w 2 X (x) (2.10) 

The kinematic variables u, 0, w, w l5 and w 2 are functions of the x coordinate only. 
Note that although the displacement assumptions in (2.1) contain seven independent 
variables, the resulting displacements given in (2.9) are in terms of only five kinematic 
variables, which are identical to those given in { 1,2} -order theory proposed by Tessler 
(1991). The terms u(x), 0(x), and w(x) are the conventional Reissner variables defined 
as weighted averages in equation (2.2). The higher-order terms for the transverse 
displacement, w { (x) and w 2 (x), account for the extension of the beam through the 
thickness. 

The quadratic transverse displacement component u z is identical to that of Tessler 
(1991). The cubic axial displacement u x has an hierarchical form such that if the higher- 
order terms Wj x and y are eliminated, the displacement field is reduced to the {1,2}- 

order theory with a linear axial displacement distribution. The further reduction of this 
theory to lower-order theories will be discussed in Section 2.2. 

The displacements given in (2.9) are the two-dimensional form of the {3,2} -order plate 
displacement field proposed by Tessler (1993b). Tessler did not formulate a {3,2}-order 
theory but proposed the use of these displacements for a hierarchical recovery of the 
{1,2} results using {3,2}-order displacement, strain, and stress approximations. This 
concept will be discussed in Section 6 as a future application of this theory. 
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2.1.2 Reduced Stress-Strain Relations 


Two sets of reduced constitutive relations have been developed for a beam with 
arbitrary material orientations: one for plane strain and the other for plane stress. For the 
development of the theory, the following notation will be utilized so that the theory is 
independent of the two-dimensional approximation employed: 


Cjj® for plane strain 

n" (k) _ 

*j ^ 

Cjj® for plane stress 


with the stress-strain relations expressed in the more general form 


^xx 

00 

C„ 

c 

'-'13 

0 " 

(k) 

r ^ 
^xx 

^zz 

> = 

C,3 

c 

'-'33 

0 

< 

Szz 

.^xz . 

g 

0 

0 

c 55 _ 

g 

y XZ; 


( 2 . 11 ) 


( 2 . 12 ) 


A complete derivation of the constitutive relations can be found in Appendix D. 


2.1.3 Strain-Displacement Relations 

Two methods are used to determine the strain-displacement relations for the beam. 

The axial strain e xx and transverse shear strain y xz are found in the usual fashion from the 
strain-displacement relations of elasticity theory. The derivation of the transverse normal 
strain departs from the conventional method. The transverse normal strain e zz is found 
by assuming a cubic stress field through the thickness for the transverse normal stress 
o zz , solving for the coefficients, and using the constitutive relations to obtain e a . The 
strain-displacement expressions are summarized here along with the of the derivation of 
the transverse normal strain. 
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The axial strain is obtained from the strain-displacement relations of linear elasticity. 


i.e., 

*^xx = U x,x = ^xO + K x0 $1 ^2 **“ K H ^3 (2*13) 

where the strain measures and curvature variables, that are functions of the x coordinate, 
and the thickness distributions are defined as 

ko • £ h]= [u„ hw Ux ] 

[ K x0> K h]=[0,x» t( W ,xx + 9,x) +W 2,xx] (2-14) 

K> ('t-'t). h[r- t) 

The transverse shear strain is obtained from the linear strain-displacement relations of 
elasticity with the additional inclusion of a shear correction factor k: 

Y“" = kYxz = k(u xz + u zx )= ky xz0 <|> xz (2.15) 

where 

[Yxz0.4>xz]=[6 + w, x! t(i-C 2 )] (2-16) 

The shear correction factor is introduced here in anticipation of correction of the shear 
stiffness of the beam. The calculation of the shear correction factor is addressed in 
Section 3. 

The transverse normal strain is not determined directly from the strain- 
displacement relations. The strain-displacement relations give rise to £ zz which is 
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continuous through the thickness, resulting in a discontinuous g zz through the thickness 
for heterogeneous beams. Analytically, the opposite is true: g^ is always continuous 
through the thickness and is discontinuous. For most cases, g Z2 is closely 
approximated by a cubic polynomial. As in Tessler (1993), the transverse normal stress 
is assumed to have a cubic expansion through the thickness: 

3 

(2.17) 

n= 0 


which leaves four coefficients G zn to be determined. Two of the coefficients are found 
using the equilibrium equation of elasticity theory, i.e., 

*,«,x + ®=,.=0 (2.1S) 


Since the transverse shear stress satisfies traction-free boundary conditions on the top 
and bottom surfaces of the beam, i.e. T xz (x,±h) = 0, the derivatives of the shear stress 
T at the top and bottom faces must vanish. To satisfy the equilibrium equation, the 
derivatives of the transverse normal stress must also vanish on the top and bottom 
surfaces: 


<J 2ZZ (x,+h)= 0 


(2.19) 


These boundary conditions reduce the transverse stress approximation to the form 

=°z0 + a z]4>5 (2.20) 

where 
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* 5 = (C-C 3 /3) 


( 2 . 21 ) 


The remaining two coefficients are found by forcing the strain field to be least- 
squares compatible with the “corrected” strain derived from the strain-displacement 
relation: 

minimize J h ( £ zz “ u z° z rr Jdz (2.22) 

with the “corrected” linear strain-displacement relation from elasticity theory given as 

(2.23) 

where k z0 and k z i are transverse correction factors. A method for obtaining these 
factors is presented in Section 3. The transverse strain measure and curvature in equation 
(2.23) are defined as 

[ezo» K zo]=[w, /h, w 2 /h 2 ] (2.24) 

At this point, must be obtained from the constitutive relations, equation (2.12), 
and is found to be 

e«=l/C«(o = -q k 3 >e x J (2.25) 

The difference between the strains defined in equations (2.25) and (2.23) is expressed as 
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A = £ (k) - u“" 


Z,Z 


^zO **" G zl 4*5 ^ 13 ) ( £ xO K xO 4*1 "*" e H 4*2 K H 4*3 ) A 

= F« ( k 

V ^ 33 


zO £ zO + ^ k zl K z0 


♦ l) 
( 2 . 26 ) 


Minimizing J A 2 dz with respect to the undetermined coefficients, g z0 and g z1 , results in 
the two equations: 


Z 0 = J_ h h A^ A dz = 0 , Z, = J_ h h A i0i] A dz = 0 (2.27) 

To simplify the equations for Z 0 and Z } , they are expressed in terms of p and r 
constants 

Z 0 = e x0 Pi + £ z oP 2 + £h p 3 +K xO p 4 +K z0 p 5 +K H p 6 +G z0 p 7 +G 2l p 8 =0 

(2.28) 

Z] = £ X o r i + £ zO r 2 + e H r 3 + K x 0 r 4 K z0 r 5 K H *6 °z0 **7 + ^ zl *8 = ^ 

where the pj and r, constants are defined in Appendix A. Equations (2.28) are solved 
simultaneously for G z0 and g zI with the results substituted into equation (2.20). 

Equation (2.25) is then simplified to yield the final expression for the transverse normal 
strain 

<4 k) = e xo¥i k) + k z0 e 2oV2 k) +e H V3 k) +KxoV4 k) + k z i K zo¥5 k) + K H V6 k) (2-29) 

where the thickness distributions \|/[ k) are defined in Appendix A. 

In contrast to the linear distribution of u™"' , equation (2.23), £ z k) has the capability of 
being discontinuous at the ply interfaces and is piece-wise cubic through the thickness. 
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This form will significantly improve the transverse strain and stress distributions through 
the thickness of the beam. 


2.1.4 Variational Principle 

The principle of virtual work is used to derive the equilibrium equations and boundary 
conditions for the beam. These equations can then be solved to obtain the five kinematic 
displacement variables which satisfy the natural boundary conditions. Neglecting body 
forces, the two-dimensional variational statement is written as 

J V (A* 5s„ + c a 84 k) + 8 y X2 )dAdx 

- j s , (q + Su 2 (x,h) )dxdy + J s _ (q“ 8u z (x 

+ I a [T x0 8u x (0,z) + T xO 8u z (0,z)]dA 

- \ A [ T xL 8 u x (L, z) + T* 8u z (L, z) ]iA = 0 


h) jdxdy 


(2.30) 


where 5 is the variational operator. A is the cross-sectional area of the beam and S + and 
S" denote the top and bottom surfaces of the beam, which, respectively, are subject to the 
normal pressure loads q + and q\ The first term is the volume integral representing the 
virtual work done by the stresses. The surface integrals denote the virtual work done by 
the external surface tractions. 

The strain-displacement relations and the displacement assumptions are substituted 
into equation (2.30) to express the virtual work principle in terms of the strains and 
kinematic variables. After integrating through the thickness and grouping terms associated 
with each virtual displacement, the variational statement may be written in terms of force 
and moment resultants as 
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J„ L [XSe*, + N z 8e^ 0 + N h 8% + M x 8k x0 + M 2 8 k z0 + M h 8k h + Q x 8y xz0 
- b(q + — q“)8w— b(q + +q”)8w 1 - b£(q + -q“)8w 2 ]dx 
+ N x0 8u(0)+ M x0 8 0(0) + M, 0 8 w u (0 ) ■ + M 2 „8 {* [e(0) + w. x (0)]+ w 2x (0)} 

+ Qx0 5w ( 0 ) + Ql0 6w i( 0 )+Q20 8w 2( 0 ) 

-N xL Su(L)- M xL 80(L)- M il Sw 1x (L)- M 2l 8{h[0(L)+ w x (L)]+ w 2x (L)} 
-QxL 8w ( L )-^lL 5w i( L )-Q2L 5w 2( L ) = 0 


The beam stress resultants N, M, and Q along with the prescribed end force and moment 
resultants N, M, and Q are defined in Appendix B. 

Since the strains and curvatures within the remaining integral are in terms of derivatives 
of the kinematic displacement variables, the stress resultants and virtual displacements are 
integrated by parts to yield 


£ [(N x ,x)8u+ (Q x ~M xx - tM H x )86 + (tM h XX -Q XX - q,)8w 
+ (N z / h + hN H xx ~ ?J 2 )8w, + (M z /h 2 +M Hxx - fq-,)5w 2 ]dx 
+ [N x „-N x (0)]Su(0)+ [M x0 -M x ( 0)]89(0) + [M l0 - hN H (0)]Sw,, x (0) 

+ [M 2O -M H (0)]6{f[e(0) + w x (o)] + W 2x (0)} 

+ |Qx»-Qx(0)]Sw(0) + [Q, 0 - h N Hx (0) ]8 w ,(0) + [Q 20 -M Hx (0)]> w 2 (0) 
-[N xL -N x (L)]8u(L) - [M xL -M x (L)] 80(L) - [M IL -hN H (L)]8w lx (L) 
-[KT 2L -Kf (l (L)]8{|[0(L) + w x ( L )]+ W 2,x(L)} 

- [QxL-Qx( L )] 8w ( L )-[QlL- hN H,x( L )] 8w i( L )-[Q2L- M H,x( L )] 8w 2( L )= 0 


where q l and q 2 are defined as 


qi(x)= b [q + (x)- q (x)] , q 2 (x)= b [q + (x) + q (x)] (2.33) 
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2.1.5 Equilibrium Equations and Boundary Conditions 


Equilibrium equations and boundary conditions are obtained from the principle of 

virtual work, equation (2.32). The expressions associated with the arbitrary kinematic 

variations must vanish independently, resulting in the following equilibrium equations: 

(Su): N x , x = 0 

(8wj): N z /h + hN Hxx — q- 2 = 0 

(Sw): fM„ >xx -Q„- q, = 0 (2.34) 

(89): Qx - M xx -f-M Hx = 0 

(8w 2 ): M z /h 2 + M Hxx - t<Ti = 0 

where the higher-order transverse equilibrium equations are associated with the 8w j and 
8w 2 variations. 

Similarly, the boundary conditions are determined by requiring that each term outside 
of the integral vanish independently. This is accomplished by either prescribing tractions 
or displacements at the ends of the beam, such that 


@x = 

= 0: 


@x = 

L: 


N x0 =N x (0) 

or 

8u(0) = 0 

n xL =n x (l) 

or 

8u(l) = 0 

M X0 =M X (0) 

or 

58(0) = 0 

m xL = m x (l) 

or 

80(L) = 0 

M I0 = hN H (O) 

or 

5w i,x(°) = ° 

M 1L =hN H (L) 

or 

8w u (L) = 0 

M 20 =M h (0) 

or 

8y(0) = 0 

-^2L == M h (L) 

or 

8y(L) = 0 

Q*o=Qx(°) 

or 

8w(0) = 0 

Qxl=Qx(l) 

or 

8 w(L) = 0 

Q,o = hN„, x (0) 

or 

8wj(0) = 0 

QlL= ^H,xW 

or 

8wj(L) = 0 

Q20 == (0) 

or 

8w 2 (0) = 0 

Q2L = -^H.x 0") 

or 

8w 2 (l) = 0 
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2.1.6 Beam Constitutive Relations 


The stress resultants which are defined in Appendix B, equation (Bl), are expanded 
using the constitutive relations of equation (2.12) and the expressions for the strains: 
equations (2.13), (2.15), and (2.29). The terms are grouped according to their associated 
quantity and expressed in terms of the beam constitutive matrix: 


fN 

X 


An 

kz(Al2 

A 13 


kziBn 

Bi 3 

0 " 


^x0 

N z 


^z0^12 

^z0^22 

k z 0^23 

^z0 ®2 1 

^zO^ zl®22 

kzoB 2 3 

0 


£ z0 

n h 


Aj 3 

kzO-A -23 

^3 3 

B3I 

^z)B 32 

B33 

0 


£ H 

•M x 

> = 

B n 

^z 0 1 

B31 

Du 

^zlDn 

D13 

0 


K x0 

M z 


k-zjBn 

^zl^zoB22 

k B 

kziDn 

k z iD 22 

^2^23 

0 


*z0 

M h 


Bl 3 

^z0 B23 

B33 

D13 

l C zlD 23 

D33 

0 


k h 

.Q*. 


0 

0 

0 

0 

0 

0 

k 2 G _ 


y xzo _ 


(2.36) 


where Ay, By, Dy, and G are the membrane, membrane-bending coupling, bending, and 
shear rigidities. These coefficients are defined in Appendix C. To complete the theory, 
the shear correction factor, k, and the transverse correction factors, k z0 and k z i, need to be 
determined for each material system investigated. The procedure for determining these 
factors is discussed in Section 3, and their numerical values are presented in Section 4. 

2.1.7 Equilibrium Equations in Terms of Displacements 

To solve the beam equilibrium equations given in (2.34), the equations must be 
expressed in terms of the five basic kinematic variables. To accomplish this, the strain 
measures and curvatures must be expressed in terms of the kinematic variables of the 
theory, i.e., 
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^xO 


1 

0 

0 

0 

0 


e z0 


0 

1 

h 

0 

0 

0 


u 

£ H 


0 

h-^r 

0 

0 

0 



K x0 

► = 

0 

0 

0 

a 

0 


w 

K z0 


0 

0 

0 

0 

1 

h 7 " 


e 

k h 


0 

0 

5 a 2 

4 9^- 

5 . a 
4^r 

a 2 

dx 2 


, W 2. 

xzO 


0 

0 

a 

1 

0 



(2.37) 


Substituting the stress resultants from equation (2.36) and the strain measures and 
curvatures from (2.37) into equations (2.34), the equilibrium equations in terms of the 
kinematic variables are given as 


(8u): A,,u xxx + kz0 h Al2 w u + hA n w Ixxx + B, ,6 xx + - z - |^’ ' 2 w 2x + 


b, 3 [r(e ,,x + w ,xxx y w 2 , x « ] = 0 


(2.38) 


(Sw, ): h {a i3 u xxx + 


k~A„ k 7l B„ 

£ W l,xx+ hA 33 W l,4x + B 3i e xxx + — 'W 2xx + 


B 33 [?( e .xxx + w ,4x X w 2 4x ] ]-f ^ {k 20 A, 2 U x + 

hk^A 23 w Ux + k z0 B 2 ,e x +^^j^-w 2 + 
k z0 B 23 [f(e,x + w„)f w 2xx ]j-- ^ 2 = 0 


^zO A 22 

h 


'Wj + 


(2.39) 


(8w): 4 ~ { B i 3 U >X xx + ■ k - Z ° h B23 W l,xx+ hB 33 W 1.4x + D 1 3^ ,xxx+ ^T^^xx + 

D 33 [t(0,xxx+ W ,4x)+ W 2.4x]j— ^ G (0 >x + W xx ) - <Tl = 0 


(2.40) 
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(§0>: B,,u*x+ + hB 31 w lxxx + D n 6 xx + 


k zl ^12 

h 2 V¥2 ’ x 


„ + 


D 13 |(e,xx + w, xxx )+w ixxx ] + f {B I3 U XX + 


k„„B 


zO 23 


W l,x + 


k„, D 


zl^23 


w,.+ 


hB,,w, + D,, 0 + = — » n 

33 l,xxx 13 ,xx j^2 2,x 

D 33[t(6,xx+ W ,xxx)+W 2 , xxx ]]- k 2 G (0 + w x )=O 


(2.41) 


(5w 2 ): B, 3 u xxx + k, °* 23 w lxx + hB 33 w 14x + D I3 6 XXX + k ^ 23 w 2 xx + 

D 33[t(0„x+ W ,4x)+ W 2.4x] + -A-{k zl B, 2 U x + Wj + 


h' 

^21^x0^22 


k 2l hB 3 2 Wi, xx +k x ,D 12 e x + 


k 2 ,D 23 [f (e, x + w xx )+ w 2 xx ]|-fq-, = 0 


w 2 + 


(2.42) 


where w 4x = w xxxx . 

Equations (2.38) - (2.42), subject to the boundary conditions given in (2.35), can now 
be solved simultaneously to determine the five kinematic variables and subsequent 
displacement, strain, and stress distributions in the beam. 

2.2 Reduction of Present Theory to Lower-Order Theories 

The hierarchical structure of the present {3,2} order theory permits a straightforward 
reduction to several lower-order theories. By eliminating the higher-order displacement 
terms w t X (x) and y(x) from equation (2.9), the displacement field reduces to the {1,2} 

form given by Tessler (1991): 
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(2.43) 


u„(x,z) = u+h£0 
u ; (x,z) = w + ^w, +(i; 2 -t)w 2 

By reducing the axial displacement to a linear approximation, the higher-order strain and 
curvature terms, £ H and k h , are eliminated from the theory. This results in the 
simplification of the equilibrium equations, boundary conditions, and stress resultants, 
given in equations (2.34), (2.35), and (2.36) respectively, such that all of the terms with a 
subscript H are eliminated. This feature will be utilized to compare the present {3,2} 
theory to the {1,2} results in Section 4. 

The {1,2} displacement theory can be further reduced to Timoshenko theory by 
neglecting the coupling between the axial and transverse stretching of the beam and by 
enforcing the inextensibility of the beam’s cross section. This is accomplished by setting 
Vj 3 = 0 and E 3 = °° . While this simplifies the equilibrium equations, the boundary 
conditions for Timoshenko theory are the same as in the {1,2} theory. The results of 
Timoshenko theory can be further reduced to those of classical beam theory by setting 
the transverse shear rigidity to be infinite, i.e., G = «> . 


3. Adjustments to Theoiy 


3.1 Shear Correction Factors 

The implementation of shear correction factors is a commonly used technique to 
correct the shear stiffness of an approximate theory. Timoshenko (1921) and Mindlin 
(1951) employed shear correction factors for the analysis of isotropic beams and plates, 
respectively. More recently, Whitney and Pagano (1970) and Reddy (1984) presented 
shear-deformable plate theories in which correction factors were applied to laminated 
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plate analysis. The aforementioned theories are based on kinematic approximations 
which tend to underestimate shear deformation, which results in an underestimation of the 
transverse displacement. The shear correction factors compensate, in a gross sense, for 
these deficiencies and bring the transverse displacement results closer to those expected 
from exact elasticity theory. 

The shear correction factor concept is introduced in the present {3, 2} -order theory to 
improve the response of laminated composite and sandwich beams. For the case of a 
homogenous beam, a shear correction factor is not needed due to the higher-order 
displacement approximations. The resulting parabolic shear stress field matches closely 
that obtained from elasticity theory. For laminated beams of arbitrary material layup, the 
shear distribution is only piece-wise smooth through the thickness. This implies that 
adjustment of the shear stiffness may be necessary and can be achieved with an 
appropriate value of a shear correction factor. 

To determine the shear correction factor within the present theory, a straightforward 
energy matching method is employed. The approach is to match the transverse shear 
energy of the present theory with that of elasticity theory for a particular “benchmark” 
problem. Since it is reasonable to expect that in the thin limit, an approximate theory 
should produce correct shear stiffness, the matching is performed for very thin beams 
(L/2h= 1000). The cylindrical bending problem, for which an exact elasticity solution is 
available, serves as our “benchmark” problem (see Section 4). Vlachoutsis (1992) and 
others have shown that, for ail practical purposes, these correction factors are 
independent of the loading and boundary conditions. The correction is therefore 
considered a function of the material system used in the analysis. This implies that once 
the correction factor is determined for a particular material system, it may be applied to 
any other analysis using that same system. Hence, it is sufficient to determine the shear 
correction factors from the cylindrical bending problem and the factors may be applied to 
other problems without loss of generality. 
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Factoring the shear correction factor from the {3,2} theory shear energy expression 
and equating it to the exact shear energy yields: 


U 


exact 

shear 


1 t t uncorrected 
k 2 shear 


(3.1) 


where is the uncorrected shear strain energy of the present theory. Both 

U uncorrected ^ 310 computed numerically. The parameter k is computed directly by 

computing the ratio of the shear strain energies. Note that the transverse correction 
factors k z0 and k z{ are also set to unity when calculating the shear correction factor as 
these factors are still unknown parameters. 

It is noted that for laminated beams, the corresponding shear correction factors are 
dependent on the material properties and layup of the plies. For each material system, a 
new shear correction factor must be computed. In Section 4.4, shear correction factors are 
computed for several material systems and layups. 


3.2 Transverse Correction Factors 

Similar to the manner by which shear correction factors improve the gross transverse 
displacement response, the transverse correction factors are used to improve the 
thickness stretch response. Mindlin and Medick (1958) showed that with each thickness 
stretch mode, i.e. constant and linear, corresponding correction factors can be specified to 
improve the transverse normal response of plates. Tessler (1995) used transverse 
correction factors for problems pertaining to the vibration of elastic plates. These factors 
were used to rectify the neglect of inertial effects in the transverse normal equations. 

While previous theories utilized correction factors to improve the dynamic response of 
plates, the present theory will extend the application of transverse normal correction 
factors to the static analysis of laminated beams. It was noted earlier that the (3, 2} -order 
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theory without the use of transverse correction factors can yield up to 50 % error for the 
following transverse normal stress boundary conditions 

(x,+h) = q + , o a (x,-h) = q" (3.2) 

These transverse normal stress boundary conditions are not explicitly enforced in the 
present theory. Because of this, the transverse stress distributions may deviate from the 
exact solution. For sandwiches with large variations in the elastic moduli through the 
thickness, these errors become even more pronounced. This error is illustrated in Figure 
3.1 for a sandwich beam in cylindrical bending in which the core properties are more than 
three orders of magnitude softer than the face sheets. The uncorrected curve (i.e., without 
the use of transverse corrections) shows the large error in the transverse normal stress on 
the top surface. 

Similar to the shear correction factors, the transverse correction factors are dependent 
on the material system, so a set of values for k z0 and k z i must be found for each layup. 
To calculate the transverse correction factors, it is convenient to introduce an alternate 
form for the beam constitutive relations, equation (2.36), relating the force and moment 
resultants to the corresponding strain measures and curvatures, i.e. 


P = [kJ[Eo][kj{^} (3.3) 

where 


A„ 

^12 

^13 

B„ 

B 12 

Dj 3 

Al2 

^22 

^23 

B 2 j 

B 22 

b 23 

A 13 

^23 
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B 3 i 

^32 

B33 


B 2 i 
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Dj 2 
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b , 2 
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®3 2 
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D 22 

D23 
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{N} t = {N x , N„ N h }, {M} t = {M x , M z , M„ }, 

= {^xO’ ^zO ’ ~ 0 » K z „ , j" 

The initial step is to obtain a solution for the case where [k z ] is an identity matrix, i.e. 
k z0 = k zl - 1 . As before, the cylindrical bending problem is analyzed in the very thin 
regime (L/2h=1000), with the appropriate shear correction factor implemented. The 
stress resultants which are obtained from this analysis are defined as {N 0 , M 0 } t . 

Substitution of these stress resultants into equation (3.3) and solving for the strain 
measures and curvatures gives rise to 

{'}-t,r[E.rM'pJ <3.4) 


such that the strain measures and curvatures are expressed in terms of the uncorrected 
computed quantities [Eq], {No}, and {M 0 } and the undetermined corrections factors k z0 
and k zl . Substituting equation (3.4) into the expressions for the axial strain (2.13) and the 
transverse normal strain (2.29) and using the stress-strain relations (2.12), the transverse 
normal stress boundary conditions can now be enforced via the transverse correction 
factors, i.e.. 


fa = (x,+h)| 

K(x,-h)j 


C+ C 3 *lje xx (k z0 ,k 2l )l 
_C' 3 C“3j|e®(k 20 ,k 2l )j 



(3.5) 


where the (+) and (-) superscripts denote the top and bottom surfaces of the beam, and 
the £ xx and strains are now functions of k z o and k zl . Solving these equations 
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simultaneously for k z o and k z i results in the numerical values for the correction factors 
that fulfill correct traction boundary conditions. By satisfying equation (3.5), the 
distribution of the transverse normal stress through the thickness closely follows that of 
elasticity theory, with correct traction conditions enforced on the surfaces. Figure 3.1 
illustrates this improvement for a sandwich beam in cylindrical bending subject to the 
normal tractions: q + = 1.0, q‘ = 0.0 (refer to Section 4, Case D, L/2h= 100 for a detailed 
description of the problem). 

This procedure is followed such that a set of correction factors is produced for each 
material system investigated. A table of transverse correction factors for the material 
systems examined in this report can be found in Section 4.4. 



Figure 3.1 Transverse Correction for a Sandwich Beam 

3.3 Integrated Interlaminar Shear Stress 

The present theory approximates the transverse shear strain y xz , equation (2.15), as a 
continuous parabolic function through the thickness. This allows accurate modeling of 
shear strain and stress through the thickness for homogeneous beams. While the parabolic 
strain satisfies the traction-free boundary conditions, the actual shear strain distribution in 
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laminated beams is often discontinuous at the ply interfaces. The parabolic 
approximation for the shear strain yields erroneous shear stress, as demonstrated in 
Figure 3.2 for an unsymmetric cross-ply graphite/epoxy laminate (refer to Section 4, Case 
C, L/2h =100). Here, the shear stress for the present higher-order theory ( HOT ) is 
obtained directly from Hooke’s law using equation (2.12), resulting in a discontinuous 
shear stress distribution. For the exact solution, the opposite is true, i.e., the shear stress 
is continuous at the ply interfaces while the shear strain is discontinuous. 




Figure 3.2 Comparison of Shear Stress and Strain for an Unsymmetric Laminate 


To improve the transverse shear stress and strain predictions, the two-dimensional 
equilibrium equation of elasticity theory 


+ T (k) = 0 

~ xz,z 


;( k) 

xx,x 


(3.6) 


is integrated to recover the transverse shear stress. The axial normal stress is 
differentiated, then integrated piece-wise through the thickness, such that the resulting 
shear stress takes the form 
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(3.7) 



a 


w 

xx,x 


dz 


The shear strain is then obtained directly from Hooke’s law: 


T 


(k) _ T (k) /p(k) 
xz — v xz ' v -"5 5 


(3.8) 


The above integration scheme is a well-established approach which generally produces 
adequate results (refer to Whitney (1972) and Reddy (1984) ). However, for 
unsymmetric and sandwich laminates, the integration approach by itself yields rather 
inaccurate shear stresses. For such laminates, the integrated shear stresses do not satisfy 
the top boundary condition even though the shear stress on the bottom surface is 
explicitly set to zero for integration. The comparison of the integrated shear stress to that 
of the exact solution leads to the conclusion that error is accumulated as a result of the 
shear stress integration through the 0° plies. The axial stretching of the midplane, a 
characteristic of unsymmetric laminates and soft-core sandwiches, introduces an error in 
these plies. The “integrated” curve in Figure 3.3 represents the integrated shear stress 
according to equation (3.7). The results are for the same material system as in Figure 3.2. 



Figure 3.3 Integration Errors Caused by Stretching of the Midplane 
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For the unsymmetric and sandwich laminates examined, the integration errors are 
distributed linearly through the 0° plies, while in the 90° plies and the core material the 
error is constant. To correct for these errors, an automated procedure is developed in 
which the total error at the top of the beam is divided into equal increments based on the 
number of 0° plies. The error is subtracted from the integrated stress curve as linear 
functions through the 0° plies, and as a constant function through the 90° plies, as 
depicted in the “Error Function” curve in Figure 3.3. Using this modified integration 
scheme, the corrected shear stresses are computed as 


T 


coir 

xz 


_ T (k) 

— L xz 


—error 

^xz 


(3.9) 


where is the integrated shear stress obtained using equation (3.7) and the error 
function is expressed analytically as 

N 

■C'=XF k (z)T k (z) (3.10) 

k=l 


where N is the total number of plies in the laminate and F k is expressed in terms of 
Heaviside functions as 

F k (z) = H(z-h k _ 1 )-H(z-h k ) (3.11) 


The shear stress error function within the ply is defined as 


* k (z) = Vi + a - 


xi 


M 


z-h 


k— 1 


V^k h-k-1 ) 


(3.12) 


where T k _j is the value of the error function at the bottom ply interface, and is initially set 
to zero ( x 0 = 0 ) at the bottom surface of the laminate. z^ n0T is the magnitude of the error 
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at the top surface of the beam, M is equal to the number of 0° plies in the laminate, 
and a is a tracer which takes the form 


jl, P = 0° 

[ 0 , p = 90° or core 


(3.13) 


where |3 is the rotation angle of the ply, as per Figure 2.2. The thickness of the ply is 
defined as t k = h k - h k _ { , where h k is the distance from the midplane to the top surface 
of the ply and h k _j is the distance from the midplane to the bottom surface of the ply. 
The variables used in equation (3.12) are depicted in Figure 3.4 for an unsymmetric 
composite beam with the lamination sequence [0/90/0/90] T . 
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Figure 3.4 Error Function and Lamination Notation 


As a result of this correction, the traction free boundary conditions are satisfied on the 
top and bottom surfaces, and the through the thickness distributions are in close 
agreement with the elasticity solution, as evident by the “corrected” curve in Figure 3.3. 
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Although not shown, the corrected shear strain is computed using equation (3.8), and 
excellent correlation to the exact solution is achieved. The shear stress distributions for 
the cases investigated in this report are in close agreement with corresponding elasticity 
solutions (refer to Section 4.5 for additional examples). Note that due to limitations in the 
exact solution. Burton and Noor (1994), only cross-ply laminates were considered for the 
cases investigated in this report. The present theory is capable of handling arbitrary 
angle-ply orientations, but the effect of these orientations on the integration scheme has 
not been investigated. 

With the exception of the homogeneous case, equation (3.9) will be used for the 
recovery of the transverse shear stress for the results presented in Section 4. It should be 
noted that the modified stresses may not exactly satisfy the original beam equilibrium 
equations. Since the modified shear stresses yield significant improvements for thick 
laminates and sandwich beams, this constraint is relaxed in order to achieve a more 
meaningful recovery of the transverse shear quantities. 
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4. Cylindrical Bending Problem 


4.1 Closed-Form Solution 

To assess the capability of the {3, 2} -order theory, a closed-form solution to the 
cylindrical bending problem is formulated and compared to the corresponding exact 
elasticity solution. Cylindrical bending is a special case of plate bending in which a sine 
load is applied to the top surface of the plate, as shown in Figure 4.1. The depth of the 
plate is infinite in the y-direction and the ends of the plate at y = 0 and y = °o are 
constrained by smooth, rigid planes such that a state of plane strain is achieved in the x-z 
plane. The plate is simply supported at the ends x = 0 and x = L such that the resulting 
deformation is uniformly cylindrical in shape. Because each cross section in the x-z plane 
deforms identically, a beam theory in a state of plane strain may be used to determine 
solutions which are identical to those obtained from plate theory. 



Pagano (1969) first published elasticity solutions for unidirectional (0°) and bi- 
directional (0°-90°) composite laminates in cylindrical bending. Burton and Noor (1994) 
developed three-dimensional elasticity solutions for laminated and sandwich rectangular 
plates. As a special case, the Burton-Noor formulation yields an exact solution for the 
cylindrical bending problem. The exact Burton-Noor solutions are used for comparison in 
the present study. 
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For the cylindrical bending problem, the loadings on the top and bottom surfaces 
of the plate are defined as 


q + (x) = q 0 sin (tc x / l) 
q"(x) = o 


(4.1) 


where q Q is the amplitude of the loading. For the solutions presented here, a unit 
amplitude will be considered. To solve the equilibrium equations derived for the present 
theory, equations (2.39) - (2.43), the kinematic displacement variables are assumed to 
have the form 


u = U cos(rcx/L), 0 = 0 cos(7tx/L) 

w = W sin(7tx/L), Wj =Wj sin(7tx/L), w 2 = W 2 sin(7tx/L) 


(4.2) 


The nature of these assumptions is such that the boundary conditions at the ends of the 
beam, equation (2.36), are all satisfied and take the form 


@x = 0: 

@ x = L: 

N x „=N x (0) 

n xL =n x (l) 

M x0 = M x (0) 

m xL = m x (l) 

M I0 = N h (0) 

m, l =n h (l) 

m 20 =m h (o) 

m 2L = m h (l) 

6w(0) = 0 

5w(L)= 0 

5w,(0) =0 

Sw/l) = 0 

8w 2 (0) = 0 

8w 2 (l) = 0 


Using the displacement assumptions in (4.2), the equilibrium equations are simplified 
such that the trigonometric functions factor out in each equation, leaving only the 
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amplitudes U, 0 , W, W ]5 and W 2 , as unknowns. Once the displacement amplitudes are 
determined, the kinematic variables are completely defined, giving rise to the strain 
measures and curvatures. The strains, stresses, and displacements are then calculated and 
plotted through the thickness for comparison to the elasticity solution. 


4.2 Material Properties 


An extended range of commonly used aerospace material systems will be investigated. 
These include a homogeneous case, such as aluminum, fiber reinforced graphite / epoxy 
composites for symmetric and unsymmetric bi-directional laminates, and sandwich beams 
composed of stiff graphite/epoxy face sheets with soft core materials. Two core materials 
are used for the sandwich cases: isotropic polyvinyl chloride foam, the most commonly 
used core material (Zenkert,1995), and titanium honeycomb, which can be modeled as an 
equivalent orthotropic material. Titanium honeycomb sandwiches are currently being 
considered by NASA for use on the High Speed Civil Transport (HSCT). The material 
properties are summarized in Table 4.1. 
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Aluminum (Al): 


Graphite / Epoxy 

(Gr/Ep): 

E l = 22.9 Msi, V LT = 0.32, G LT - 0.86 Msi, 

E T — 1.39 Msi V TT = 0.49, G XT — 0.468 Msi 

Polyvinyl Chloride 
(PVC): 

E = 15.08 ksi, V = 0.3, G = 5.80 ksi 

Titanium Honeycomb 
(Ti): 

B 1 = 62.36 psi, V 13 = 5.6x1 O' 5 , G 13 = 75.1x1 0 3 

psi, 

E 2 = 41.27 psi, V 23 = 3.7x1 O' 5 , G 23 = 56.7x1 0 3 

psi, 

E 3 = 345x1 0 3 V 12 = 1.23, G 12 = 1140 psi 

psi, 


L = Longitudinal Direction (fiber), T = Transverse 
1,2, & 3 = principle material directions 


Table 4.1 Material Properties 


4.3 Test Case Definition 

Six test cases are investigated using the material systems given in Table 4.2. While 
isotropic materials have identical material properties regardless of the orientation in the 
coordinate frame, orthotropic materials, like graphite/epoxy, are commonly oriented in 
various directions in the x-y plane to achieve the desired stiffness characteristics. The 
layup column in Table 4.2 gives the order of lamination, beginning with the bottom ply at 
z = -h , and the angular orientation of the longitudinal fiber direction for the Gr/Ep plies 
with respect to the x axis. The S subscript denotes symmetry with respect to the 
midplane of the beam. For unsymmetric laminates, the total layup must be given and is 
denoted by the T subscript. The unsymmetric cases will test the coupling effect between 
the stretching and bending modes of the beam. Because of the theoretical restrictions on 
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the exact elasticity cylindrical bending solution, only bi-directional 0° and 90° orientations 
are used in the analysis. 


Case 

Material 

Layup 

A: 

A1 

Homogeneous 

B : 

Gr/Ep 

[ 0 2 /90 2 /0 2 /90 2 3s 

C: 

Gr/Ep 

[ 0 4 /90 4 /0 4 /90 4 ] T 


Gr/Ep - PVC 

[ 0 4 /90 2 /0 4 /90 2 /0 4 / PVC Core ] s 


Gr/Ep - Ti 

[ 0 4 /9 0 2 /0 4 /9 0 2 /0 4 / Ti Core ] s 


Gr/Ep - PVC 

[ 0 4 /90 2 /0 4 /90 2 /0 4 / PVC Core 
/90 4 /0 2 /90 4 /0 2 /90 4 ] T 


Table 4.2 Definition of Test Cases and Corresponding Ply Layup 


The ply thickness, t k , for each graphite/epoxy lamina is set equal to 0.00625 in. This 
is a reasonable assumption for the actual thickness used in the manufacturing of 
composite laminates. The total thickness of the beam for the isotropic Case A and the 
sandwich Cases D, E, and F is 1 .0 inch. The cores for the sandwich cases comprise 80% 
of the total thickness of the beam. For the composite laminate. Cases B and C, the total 
thickness is equal to 0.1 inch. The length of the beam will be varied to achieve the desired 
span to thickness ratio. 

4.4 Correction Factors 

The analytic solutions of the {3,2} theory rely on the use of appropriate shear and 
transverse correction factors. Following the procedure in Section 3, the correction factors 
are calculated for very thin beams (span to thickness ratio L/2h = 1000) and are 
summarized in Table 4.3. 
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Case 

k 2 

k Z 0 

kzi 

A 

1.0 

1.0 

1 .0 

B 

0.91056 

1.0 

0.96494 

C 



1.01975 

D 

0.30666 

1.24326 

1.59569 

E 

0.85883 

1.14298 

1.26431 

F 

0.33811 


1.59787 


Table 4.3 {3,2} Theory Correction Factors 


Note that the unsymmetric Cases, C and F, require a significant amount of shear 
correction as compared to the symmetric cases. This is because non-symmetry causes a 
coupling effect between the stretching and bending modes of the beam, resulting in a 
significant stretching of the midplane. 

To afford a numerical comparison to the {1,2} -order theory, the corresponding 
correction factors for that theory are also computed and summarized in Table 4.4. 


Case 

k 2 

k z o 

kzi 

A 

1.0 

1 .0 

1.0 

B 

0.97793 

1.0 

1.0 

C 

0.73262 

1.0 

1.0 

D 

0.37301 

1.24326 

1.76405 

E 

1.04737 

1.14298 

1.39621 

F 

0.41210 

1 .24326 

1.76405 


Table 4.4 {1,2} Theory Correction Factors 


4.5 Results 

Six test cases are investigated over a full range of span to thickness ratios. For each 
case defined in Table 4.2, thin, moderately thick, and very thick beams with the span to 
thickness ratios L/2h ={ 100, 10, 4 } are analyzed. A unit beam width is assumed. The 
displacement, stress, and strain distributions through the thickness of the beam are 
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calculated at the locations where the quantities are maximum. These locations are 
depicted in Figure 4.2. 


* L * 


@ x=0 

@ x=L/2 

Bx ’ ^xz > Yxz 

> Fjk , Fzz , 5 (j 22 


~T 

2h 

Jl 


Figure 4.2 Locations for Displacement, Strain, and Stress Computations 


Because the axial strain £ xx is simply the derivative of the axial displacement u x , their 
maximum values are proportional according to the equation £ xx (L/2) = - u x (0) . The 
shear strain y xz is not shown in the graphical results since it is directly proportional to 
the shear stress r xz , as per equation (3.8). While the strains are dimensionless, it is also 
desirable to have the transverse displacement and stresses expressed in terms of non- 
dimensional quantities for a more general comparison of results. The transverse 
displacements are normalized with respect to the exact deflection at the midplane 
u z (£ = 0) exact . This normalization will also emphasize the nonlinear distribution of the 
transverse displacement for thick beams. The stresses are normalized with respect to the 
applied normal traction q, which at mid-span of the beam is equal to unity. 

For each span to thickness ratio, six plots are presented to compare the present higher- 
order theory to the exact solution. For the smallest aspect ratio, L/2h = 4, the present 
{3,2} theory is also compared to the {1,2} theory. As the span to thickness ratio 
approaches the thin regime, the {1,2} theory results are expected to approach those 
obtained by the present theory, so a {1,2} comparison is not shown for span to thickness 
ratios greater than 4. For the test cases, the following notation is used for graphical 
comparison of the theories: Exact stands for the exact elasticity solution, HOT denotes 
the present {3, 2} -order theory, and { 1 , 2 } represents the results of the {1,2} -order 
theory. 
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4.5.1 Homogeneous Beams 


Case A assesses the accuracy of the theory for analysis of homogeneous materials. 

As demonstrated in Figures 4.3 - 4.5, the displacements, strains, and stresses, for both the 
{3,2} and {1,2} theories, produce excellent correlation with exact elasticity solutions. 

This is accomplished without the use of shear or transverse correction factors. Notice 
that as the span to thickness ratio decreases, the transverse deflection takes on a parabolic 
distribution through the thickness. For very thick beams (L/2h= 4), this deflection varies 
up to 3% of the nominal midplane deflection. The axial displacement, strain, and stress, 
as well as the transverse strain all take on linear distributions through the thickness 
regardless of the span to thickness ratio. At this point, the present {3,2} displacement 
theory does not provide any benefit over the {1,2} theory, but the results shown here are 
presented to demonstrate the accuracy of these theories for commonly used homogeneous 
materials. 

The results in Figures 4.3 - 4.5 are expected since the exact normal strain and shear 
strain distributions are closely approximated by linear and parabolic distributions, 
respectively. These distributions are easily captured with only a linear axial displacement 
assumption and a quadratic transverse displacement. Note that both {3,2} and {1,2} 
theories have a unique advantage over previous higher-order shear deformation theories: 
they satisfy the traction-free boundary conditions without the need to integrate the 
equilibrium equations of elasticity. Many of the higher-order deformation theories do not 
satisfy the traction-free boundary conditions and cannot produce accurate shear stress 
results without integrating the equilibrium equations of elasticity theory. Though not 
shown, results of similar accuracy can be obtained for homogeneous orthotropic beams 
without the use of correction factors, i.e., k 2 = k z o = k z i = 1.0. 
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Figure 4.4 Case A: Aluminum, L/2h = 10 
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HOT & {1,2} 



Figure 4.5 Case A: Aluminum, L/2h 








4.5.2 Laminated Composite Beams 


Figures 4.6 - 4.8 depict through the thickness distributions of displacements, strains, 
and stresses for very thin (L/2h= 100), moderately thick (L/2h= 10), and very thick (L/2h 
= 4) symmetrically laminated composite beams. The material system is that identified in 
Table 4.2 as Case B. Similar distributions for an unsymmetric laminate, Case C in Table 
4.2, are shown in Figures 4.9 - 4.1 1. 

Figure 4.6 depicts displacement, strain, and stress results for a symmetrically 
laminated, thin beam. The {3,2} results are in excellent agreement with the exact 
solutions. While the transverse normal stress varies slightly from the exact curve through 
the thickness, the midplane stress and surface tractions are all accurately recovered. In 
Figure 4.7, the thickness distributions for the moderately thick case remain in close 
agreement with the elasticity solution. Figure 4.8 compares the results of the {1,2}, 

{3,2}, and exact solutions in the very thick regime. At the midplane, the transverse 
displacement is within 1% of the exact solution for both approximations, whereas at the 
top surface the difference is about 6%. The {3,2} cubic approximation of the axial strain, 
£ xx , produces accurate maximum values at the top and bottom surfaces, whereas the 
{1,2} linear distribution underestimates these maximum values. This results a 20% 
underestimation for G xx when the {1,2} -order theory is used. This discrepancy is 
especially significant when designing a laminate according to the maximum stress criterion. 

For the case of an unsymmetric laminate, the midplane exhibits straining due to the 
axial and transverse strains and the coupling between the two stretching modes. In Figure 

4.9, a tensile e xx and a compressive are observed at the midplane. In Figures 4.9 and 

4.10, close agreement is shown for the through the thickness distributions, with the 
exception of the transverse normal stress a zz . The midplane stretching results in a G zz 
thickness variation that is not captured by the higher-order theory. In Figure 4.1 1, the 
{3,2} theory shows improvements over the {1,2} theory for very thick beams, although 
the polynomial distributions are unable to capture the piece-wise smooth, exact solution. 
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Figure 4.6 Case B: Gr/Ep Symmetric Laminate, L/2h = 100 
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Figure 4.7 Case B: Gr/Ep Symmetric Laminate, L/2h = 10 
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Figure 4.8 Case B : Gr/Ep Symmetric Laminate, L/2h = 4 
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Figure 4.9 Case C: Gr/Ep Unsymmetric Laminate, L/2h =100 
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Figure 4.10 Case C: Gr/Ep Unsymmetric Laminate, L/2h =10 
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Figure 4.11 Case C: Gr/Ep Unsymmetric Laminate, L/2h = 4 
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4.5.3 Sandwich Beams 


Cases D, E, and F address the analysis of symmetric and unsymmetric sandwich 
composite beams. Sandwich beams present a unique challenge to the theory owing to the 
drastic change in the material properties through the thickness. The face sheets of a 
sandwich are stiff in the axial direction to carry the axial stress while the core material 
may be considerably less stiff. This is possible because the magnitude of the transverse 
shear load is usually much less than the axial load. As a result, the core can be made of a 
lightweight, compliant material, reducing the weight of the structure. Typical thickness 
dimensions for the components of a sandwich beam are illustrated in Figure 4.12 



Figure 4.12 Sandwich Composite Beam 


Figure 4.13 shows the results for Case D in the thin regime. Note how the axial stress 
G xx is carried by the Gr/Ep face sheets on the top and bottom faces of the beam while the 
transverse shear stress T xz is almost exclusively carried by the PVC core. For the thin 
sandwich, the {3,2} -order theory almost exactly matches the exact elasticity solution, 
while for the thicker case, shown in Figure 4.14, the limitations in the theory are more 
evident. The transverse displacement prediction is about 2% in error. The stresses and 
strain on the top and bottom faces, where these quantities are usually the largest, are 
accurately predicted by the theory. For the thick regime. Figure 4.15, the transverse 
displacement is overestimated by both the {3,2} and {1,2} theories, with errors ranging 
from 12% to 37%. Notice that the linear approximation for the axial displacement in the 
{1,2} theory underestimates the axial strain, resulting in a significant underestimation of 
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the axial stress. The {3,2}-order theory captures G xx at the top and bottom surfaces 
adequately, while the same stress for the {1,2} theory is 85% in error. The quality of 
these results suggests that the span to thickness ratio of L/2h = 4 may constitute a 
practical limit of application of this theory for sandwich beams. 

In Case E, the sandwich core is made of titanium honeycomb. Because of the 
orthotropic stiffness properties of the honeycomb, the axial stiffness is very low 
compared to the relatively high transverse stiffness. The sandwich core acts as a shear 
web, resisting the transverse stresses only. The titanium honeycomb serves as an efficient 
core material, providing a large transverse stiffness with a minimal weight penalty. Figures 
4.16 through 4.1 8 depict the results for the sandwich beam with a titanium honeycomb 
core. Note that the nonlinear axial effects are present until the beam approaches the thick 
regime (L/2h = 4). Comparing these results to those for the PVC core shows that the axial 
displacement and strain are already nonlinear at L/2h =10. The axial quantities become 
increasingly nonlinear as the differences between the transverse stiffness through the 
thickness of the laminate increase. This is evident with the larger discrepancy between 
the {3,2} and {1,2} results for the sandwich beam with a PVC core. As a result, the 
higher-order terms in the {3,2} theory do not provide as much benefit over the {1,2} 
theory for Case E with a titanium honeycomb core. 

In Case F, an unsymmetric sandwich beam with a PVC core is analyzed. This 
particular sandwich presents the greatest challenge to the beam theory because of the 
axial-bending coupling and the low stiffness of the PVC core. Figures 4.19 through 4.21 
depict accurate correlations with the exact solution. In Figure 4.21, the benefits of the 
higher-order axial displacement over the linear axial displacement are evident. While the 
deflection preducted by the {3,2} theory differs from the exact solution by about 8%, the 
{1,2} theory prediction is in error by 22%. The axial strain exhibits a highly nonlinear 
distribution through the thickness, with the boundary conditions accurately modeled with 
the {3,2} theory. Although not particularly evident from Figure 4.21, the axial 
compressive stress for the {1,2} theory is underestimated on the bottom surface by 83%. 
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Figure 4.13 Case D: Gr/Ep-PVC Symmetric Sandwich, L/2h = 100 
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Figure 4.14 Case D: Gr/Ep-PVC Symmetric Sandwich, L/2h = 10 
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Figure 4.15 Case D: Gr/Ep-PVC Symmetric Sandwich, L/2h = 4 
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Figure 4.16 Case E: Gr/Ep-Ti Symmetric Sandwich, L/2h = 100 
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Figure 4.17 Case E: Gr/Ep-Ti Symmetric Sandwich, L/2h = 10 
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Figure 4.19 Case F: Gr/Ep-PVC Unsymmetric Sandwich , L/2h = 100 
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Figure 4.20 Case F: Gr/Ep-PVC Unsymmetric Sandwich , L/2h =10 
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Figure 4.21 Case F: Gr/Ep-PVC Unsymmetric Sandwich , L/2h = 4 
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5. Guidelines for Finite Element Approximations 


To enable the use of the present {3,2} theory on a wider range of applications, finite 
element approximations of the theory need to be undertaken. Tessler (1991) developed a 
computationally attractive two-node beam finite element based on the { 1,2} -order theory. 
This element has the same number of degrees of freedom (dof) as a comparable 
Timoshenko beam element, yet is capable of recovering higher-order displacement and 
stress distributions through the thickness. While Timoshenko theory only accounts for 
transverse shear deformation, the {1,2} -theory element captures both the transverse shear 
and transverse normal deformations without any additional computational cost. The 
higher-order displacements w, and w 2 are statically condensed out at the element 
equilibrium level, leaving only u, 0 , and w as kinematic degrees of freedom. Thus, the 
resulting element has only 6 dof, 3 dof for each node. 

It is noted that thick composite and sandwich structures usually exhibit nonlinear axial 
displacement, strain, and stress distributions. For such structures, application of finite 
elements derived from the {1,2} theory may result in somewhat inadequate predictions. 
Instead, finite elements based upon the {3,2} theory can be used to model regions where 
the beam is particularly thick. This will allow the higher-order stress effects to be 
accurately modeled while keeping the simplicity of the beam element. Also, the basic 
kinematic variables u, 0, and w in the {3,2} theory are have the same physical meaning 
as in the {1,2} theory, ensuring proper compatibility between the elements. 

Examining the highest order derivatives of the kinematic variables in the variational 
statement (2.31), necessary continuity requirements for the kinematic variables can be 
established (Hughes, 1987). While the shape functions within the finite element are 
always required to be continuous, the shape functions at the element interfaces, i.e. nodes, 
must satisfy minimum continuity requirements. For the kinematic variables which appear 
as first derivatives in the variational statement, i.e. u x and 0 x , the shape functions must 
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be C° continuous, meaning they must be continuous at the element nodes. For kinematic 

variables which appear as second derivatives in the variational statement, 

i.e. w xx , w lxx , and w 2 xx , the shape functions must be C 1 continuous. This implies that 

both the function itself and its first derivative must be continuous. To accomplish this, 
degrees of freedom at each node must be associated with the variables: 

u, 0, w, Wj, and w 2 
w ,*, w i,x - and w 2x 

The {1,2} theory enables the Wj and w 2 higher-order variables to be condensed out at 
the equilibrium equation level, thus defining w, and w 2 in terms of u and 0, respectively. 
Since only u, 0 , and w need to be computed at each node, this allows for a simple, two- 
node, 6 dof element. A two-node element based on the present {3,2} theory will have 16 
dof, 8 dof at each node, and the number of degrees of freedom cannot be reduced. 

Considering a relatively large number of dof for a {3, 2} -order beam finite element, a 
global-local computational strategy can be invoked. For example, the {1,2} beam 
elements would be sufficient to model the majority of beams, whereas the {3,2} elements 
would be used locally in thick regions of a beam where nonlinear effects are present. An 
example of this global-local modeling strategy would be a tapered laminated beam where 
the thickness increases linearly along the span. Figure 5.1 demonstrates the modeling of 
such a beam using both {1,2} and {3,2} elements. 

It is noted that, as in the {3,2} analytic theory, the transverse shear and transverse 
normal correction factors are inherently present in the element stiffness matrix. Since the 
correction factors are only dependent on the material system, one set of factors is needed 
for the material layup used in the analysis. As previously discussed, the cylindrical 
bending problem provides a convenient means of computing the correction factors. The 
associated calculations are easily automated and computationally very efficient. 


64 




{1,2} Beam Elements 


• 3 dof per node 


Compatible Element Interface 



{3,2} Beam Elements 




• 8 dof per node 


Figure 5.1 Finite Element Model of a Tapered Beam 
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6. Conclusions and Recommendations 


A higher-order bending theory for laminated composite and sandwich beams has been 
successfully developed. Using a special cubic distribution for the axial displacement, the 
{1,2} theory has been extended to include higher-order axial effects without introducing 
additional kinematic variables. A special quadratic form of the transverse displacement is 
assumed such that the traction-free shear stress boundary conditions on the surface of the 
beam are satisfied. An independent expansion is also assumed for the transverse normal 
stress, notably improving the transverse strain and stress predictive capability of the 
theory. Appropriate shear correction factors have been used to adjust the shear stiffness 
based on energy considerations, thereby improving the transverse displacement response. 
A unique implementation of transverse correction factors has led to significant 
improvements in the transverse stretch response for laminated composite and sandwich 
beams. Accurate shear stresses for a wide range of laminates, including the challenging 
unsymmetric laminates and sandwich composites, have been obtained using an original 
integration scheme, which corrects for errors introduced due stretching of the midplane. 

A closed-form solution to the cylindrical bending problem has been derived, 
demonstrating excellent correlation to elasticity solutions for a wide range of beam aspect 
ratios and commonly used material systems. Results show that the axial response for 
thick sandwich beams has been improved over the { 1 ,2} theory. To enable the use of the 
present {3,2} theory on a wider range of applications, guidelines for finite element 
approximations have been presented. 

There are several possibilities for future work based on this theory. First, a 
straightforward extension of the present theory to plate theory should be possible, 
thereby enabling the analysis of laminated composite and sandwich plates. Second, a 
finite element formulation could be implemented to assess the accuracy of {3,2} -order 
beam elements. Although the kinematics of the {3,2} theory would require more nodal 
degrees of freedom than traditional finite elements, the higher-order beam elements may 
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provide improved accuracy previously attained only by complex three-dimensional finite 
elements. 

Since the analytical theory and corresponding finite element approximation are more 
complex than the those of the {1,2} theory, it may be possible to employ a more efficient 
computational strategy while retaining the accuracy of the {3,2} theory. One possibility 
would be to determine the five kinematic variables from the {1,2} theory and then 
substitute these variables into the {3,2} equations to obtain displacement, strain, and 
stress distributions. Since the kinematic variables have the same meaning in each theory, 
these variables should be comparable in magnitude. This concept has been proposed by 
Tessler (1993b) with promising results for thick laminates, demonstrating significant 
improvement in the predictive capability of the {1,2} theory using the {3, 2} -order 
recovery. 
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Appendix A: Transverse Normal Strain Coefficients 


The shape functions for the transverse normal strain defined in equation (2.29) are defined 
as 

v! k) = S<*> (p, -R, 0 5 -c< 3 k > ) = S< k > (p 4 - r 4 <^ 5 -C,f <(>, ) 

= S< k > (p 2 - R 2 <t> 5 ) v' 5 l ‘ = SS> (p 5 -R 5 <t) 5 ) (Al) 

Vf = S< k 3 > (P 3 - R 3 <t> 5 - C< k > <|) 2 ) V < k) = S« k » (P 6 - R 6 <|> 5 - C« «t>, ) 
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Appendix B: Stress Resultants and Prescribed Tractions 


Beam Stress Resultants: 
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Prescribed Stress Resultants at Ends of Beam (@x=0 & x=L ) : 
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Appendix C: Stiffness Coefficients for Beam Constitutive Matrix 


Ay Stiffness Coefficients: 
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Appendix D: Derivation of Constitutive Relations 

The stress-strain relations for a two-dimensional stress state are developed for plane 
strain and plane stress in the x-z plane (refer to coordinate system in Figure 2.1). Both 
constitutive relations are derived from the general three-dimensional state of stress and 
strain for an orthotropic ply. The constitutive relations for the k th ply of a laminated 
orthotropic solid are defined as 
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which in condensed form are 

{*.r- i.rt>.r 


where the (k) superscript implies that the quantity is a function of the k^ ply; [S m ] is 
the conventional orthotropic compliance matrix, and the m subscript identifies the lamina 
axes in the material coordinate system. The vectors { o m } and { £ m } contain stresses and 

strains defined in the material reference frame. 

The orthotropic plies of a laminated composite or sandwich may be of arbitrary 
orientation in the x-y plane of the beam. Since the constitutive relations are given with 
respect to the material coordinate system, the properties must be transformed into the 
beam coordinate system. The following matrices define the stress and strain 
transformations : 
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c = cosp , s = sin|3 


where |3 denotes the angle between the principal material direction of the ply and the x 
coordinate, as depicted in Figure D.l. 



Figure D. 1 Beam and Material Coordinate Frames 

The stress and strain quantities for the material reference frame are expressed in terms 
of the global (x-y coordinate) stresses and strains using the transformation matrices 
defined in (D2): 

KTMtJW" • fc.f-WfeF < m > 

where the g subscript identifies quantities defined in the global coordinate system. The 
equations (D3) are substituted into (Dl) to yield 
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(D4) 


krt.r- [sj k, [T 0 ]^ g f 


Pre multiplication of both sides of the equation by the inverse of the strain transformation 
matrix results in the form 

« k, = &nj <P5) 

with the transformed compliance matrix defined as 


Expressing equation (D5) in expanded form yields the transformed three-dimensional 
strain-stress relations in terms of global quantities: 
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Although the same notation is used, the Sy coefficients shown in (D7) are not the same as 
those in (Dl). The coefficients from equation (D7) are derived from the transformed 
compliance matrix defined in (D6) and are not expanded individually for this discussion. 
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Pre multiplication of equation (D7) by y = [C g ] results in the three- 
dimensional stress-strain relations: 



where the stiffness matrix is defined as 

[cj 1 - M' [C.fM <D9> 

Plane Strain 

For the case of plane strain in the x-z plane, the strain components e y , y y2 , and y xy 

are assumed to be negligible. By setting these strains equal to zero and solving equations 
(D7) simultaneously, the reduced constitutive relations for plane strain become 
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The plane strain constitutive relations can be inverted to yield the stress-strain relations: 
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Plane Stress 

For a state of plane stress in the x-z plane, the stress components o y , x yz , andx xy in 

equation (D8) are set equal to zero. Solving the equations in terms of the stress 
components associated with x and z yields the two-dimensional constitutive relations for 
plane stress: 
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